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Computer  moleculiar  dynamics  simulation  studies  of  grain-boundary 
structures.  II.  Migration,  sliding,  and  annihilation  in  a  two-dimensional  solid 
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02139 

(Received  22  December  1981;  accepted  for  publication  7  April  1982) 

The  dynamics  of  a  high-angle  tilt  boundary  in  a  triangular  lattice  with  a  truncated  Lennard-Jones 
potential  has  been  simulated  by  the  method  of  computer  molecular  dynamics.  Boundary 
migration  was  observed  at  various  temperatures  above  a  threshold,  and  was  interpreted  in  terms 
of  small  steps  involving  coupled  sliding  and  migration.  In  some  cases,  the  simulation  proceeded  to 
the  annihilation  of  a  pair  of  boundaries  with  the  attendant  release  of  the  boundary  free  energy  and 
delocalization  of  the  boundary  free  volume. 

PACS  numbers:  61.70.Ng,  66.30.Lw,  62.20.Fe 


I.  INTRODUCTION 

Despite  the  important  role  of  grain  boundaries  in  mate¬ 
rials  properties,  our  knowledge  of  how  boundaries  actually 
move  at  the  microscopic  level  is  limited. 1  Much  of  the  diffi¬ 
culty  is  due  to  the  lack  of  a  suitable  means  of  observing  the 
dynamical  processes  with  sufficient  spatial  and  time  resolu¬ 
tion.  However,  through  computer  simulation,  one  can  ob¬ 
tain  detailed  atomistic  information  about  model  systems, 
and  depending  on  the  choice  of  interatomic  potential  func¬ 
tion,  it  is  possible  to  simulate,  to  a  quite  reasonable  level  of 
accuracy,  structural  and  dynamic  properties  of  real  materi¬ 
als.2-*  Both  static7  and  dynamic8-11  simulations  of  grain 
boundaries  have  been  carried  out.  l.i  the  latter  studies,  one 
follows  the  atomic  trajectories  by  numerically  integrating 
the  Newton’s  equations  of  motion  over  a  succession  of  small 
time  intervals  for  a  system  consisting  of  several  tens  to  se¬ 
veral  hundreds  of  particles. 

In  this  paper,  we  present  results  on  computer  molecular 
dynamics  simulations  of  a  two-dimensional  solid  containing 
two  high-angle  tilt  boundaries.  The  system  was  a  cell  of  56 
atoms  interacting  via  a  Lennard-Jones  6-12  potential,  and 
periodic  border  conditions  were  employed  to  eliminate  sur¬ 
face  effects.  Initially,  the  system  was  relaxed  at  zero  tem¬ 
perature,  then  it  was  brought  up  to  equilibrium  at  various 
temperatures  and  the  behavior  of  the  grain  boundaries  moni¬ 
tored  by  following  the  atomic  positions  as  a  function  of  time. 

The  principal  result  of  the  study  is  the  observation  of 
the  detailed  atomic  motions  involved  in  grain  boundary  mi¬ 
gration.  At  sufficiently  high  temperature  the  two  boundaries 
in  the  simulation  cell  both  underwent  migration  and  even¬ 
tually  annihilated  to  form  a  single  crystal.  The  migration 
was  found  to  be  thermally  activated,  since  below  a  certain 
threshold  temperature,  no  migration  occurred  in  the  time 
interval  of  simulation. 

The  present  results  therefore  indicate  that  atomic  mo¬ 
tions  in  a  grain  boundary  are  highly  cooperative  and  depend 
to  a  considerable  extent  on  the  boundary  structure.  More¬ 
over,  they  show  that  migration  involves  a  local  shear  at  the 
boundary,  thus  implying  that  grain  boundary  migration  in 


real  solids  occurs  by  a  grain  boundary  dislocation  mecha¬ 
nism.  We  will  also  see  below  that  the  concept12  of  the  DSC 
(Displacement  Shift  Complete)  lattice  is  particularly  useful 
in  interpreting  the  computer  results. 

II.  CONSTRUCTION  OF  SIMULATION  CELL 

The  simulation  cell  contains  all  the  particles  whose  tra¬ 
jectories  are  to  be  calculated  explicitly  using  the  equations  of 
motion.  Periodic  border  conditions  have  the  effect  of  sur¬ 
rounding  the  simulation  cell  with  image  cells  generated  by 
translating  the  simulation  cell  along  the  x  an&y  directions  by 
distances  equal  to  the  cell  dimensions  in  the  two  directions. 
They  eliminate  explicit  surface  effects  by  ensuring  that  parti¬ 
cles  in  the  simulation  cell  close  to  the  cell  border  will  see  the 
same  environment  as  if  they  were  in  the  cell  center.  What¬ 
ever  defect  structure  that  is  present  in  the  simulation  cell  is 
reproduced  in  the  image  cells,  and  the  simulation  is  that  of 
an  infinite  array  of  defects.  Normally  one  needs  to  consider 
only  the  immediate  neighboring  image  cells  when  consider¬ 
ing  forces  on  atoms  in  the  simulation  cells,  but  when  atoms 
move  more  than  a  small  fraction  of  the  cell  dimension,  it  is 
necessary  to  include  more  distant  image  cells.1* 

The  simulation  cell  is  constructed  by  considering  a  bi- 
crystal  of  two  identical  hexagonal  lattices  which  are  rotated 
with  respect  to  one  another  as  shown  in  regions  A  and  B  of 
Fig.  1(a).  The  rotation,  in  this  case  38.21*  between  close 
packed  directions,  is  chosen  to  be  a  coincidence  misorienta- 
tion;  that  is,  one  which  brings  a  certain  fraction,  1/2,  of  the 
lattice  sites  of  the  two  rotated  lattices  into  geometric  coinci¬ 
dence  (2  =  7).  This  can  be  seen  in  region  C  where  the  two 
crystals  interpenetrate.  The  coincidence  sites  form  a  super- 
lattice,  the  coincidence  site  lattice  ICSL).M  A  coincidence 
boundary  was  chosen  for  study  because  the  boundaries  are 
fundamental  to  current  geometric  models  for  grain  bound¬ 
ary  structure. 121,16  Further,  such  boundaries  have  a  period¬ 
ic  structure,  the  periodicity  of  the  CSL  giving  the  periodicity 
of  possible  grain  boundaries  at  the  coincidence  misorienta- 
tion,  consistent  with  the  use  of  periodic  borders  on  cell  faces 
normal  to  the  boundary  plane. 
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Figure  1(a)  also  shows  in  region  D  a  fine  rectangular 
grid  delineating  the  { 541 )  and  j  321  j  planes  in  each  crystal. 
Half  the  points  in  this  grid  also  represent  points  in  the  DSC 
lattice. 1 2  To  understand  its  significance,  one  should  note  that 
a  general  translational  displacement  of  one  of  the  rotated 
lattices  with  respect  to  the  other  will  disrupt  the  coincidence 
pattern  and  the  coincidence  sites.  However,  when  this  dis¬ 
placement  is  a  DSC  lattice  vector,  the  coincidence  pattern 
will  be  restored,  usually  with  a  shift  in  location  of  the  bound¬ 
ary.  In  the  following,  the  DSC  lattice  is  repeatedly  used  as  a 
metric  in  constructing  the  computational  cell  and  in  describ¬ 
ing  atom  displacements  as  Fig.  1(a)  suggests.  More  impor¬ 
tantly,  the  lattice  plays  a  fundamental  role  in  the  interpreta¬ 
tion  of  the  migration  behavior  observed  in  the  simulation. 

The  DSC  lattice  for  the  2  =  7  CSL  is  also  hexagonal 
with  a  nearest-neighbor  distance  1/7  that  of  the  CSL.  Again 
it  is  more  convenient  to  use  a  centered  rectangular  cell.  The 
unit  vectors  of  that  cell  are  [100]  DSC  =  1/7  [100] 
CSL  =  (a/14)  <321), and [010]  DSC  =  1/7  [010] CSL  =  (a/ 
14)  (541).  The  detailed  relation  between  the  fine  grid  in  Fig. 
1(a)  and  the  hexagonal  and  rectangular  DSC  unit  cells  are 
shown  in  Figs.  1(b)  and  1(c). 

Our  molecular  dynamics  simulation  was  initialized  us¬ 
ing  a  statically  relaxed  cell,  with  dimensions  XL  =  7  [100] 
DSC  =  CSL  in  the  x  direction  and  YL  =  28.607  [010] 
DCS  =  4.087  CSL  in  the  y  direction.  Dimension  YL  was 
determined  by  adjusting  the  free  volume  at  the  boundary,  as 
will  be  discussed  later.  The  56-atom  simulation  cell  contains 
two  grain  boundaries,  the  second  boundary  being  a  result  of 
the  periodic  border  conditions  on  they  border;  for  example. 
Fig.  1(a)  shows  the  bottom  half  of  the  kite  structure  in  the 
bottom  border  of  the  cell.  The  two  boundaries  are  crystallo- 
graphically  identical  but  of  opposite  rotational  sense.  The 
atoms  in  the  two  boundaries  are  out  of  the  range  of  direct 
interaction  through  the  pair  potential;  however,  the  boun¬ 
daries  affect  each  other  slightly  through  the  long-range  com¬ 
ponents  of  their  elastic  strain  field.  Because  of  the  opposite 
rotational  sense,  one  may  expect  that  if  they  are  moved  close 
enough  to  each  other  they  will  attract  since  their  mutual 
annihilation  would  reduce  the  system  potential  energy. 
From  the  present  results  it  appears  that  the  attraction  per¬ 
sists  even  out  to  their  initial  separation  distance.  Since  the 
boundaries  initially  form  part  of  an  equally  spaced  infinite 
array,  there  must  be  zero  force  on  any  individual  boundary 
from  its  neighboring  boundaries  until  the  symmetry  is 
broken  by  the  motion  of  a  boundary  with  respect  to  its  near¬ 
est  neighbors.  In  the  simulation  the  two  boundaries  are  dy¬ 
namically  independent,  in  the  sense  that  atom  motions  at  the 
two  boundaries  are  uncorrelated  except  for  long-time  aver¬ 
age  behavior  reflecting  the  elastic  interaction. 

The  particular  YL  value  above  was  reached  after  a  stat¬ 
ic  relaxation  study  to  obtain  a  0'K  equilibrium  structure 
with  the  proper  free  volume  Vf  at  the  boundary. 17  Here  V,  is 
the  volume  occupied  by  atoms  in  a  box  containing  a  grain 
boundary  in  excess  over  the  volume  occupied  by  the  same 
number  of  atoms  in  a  perfect  crystal,  the  extent  of  the  box  in 
the  direction  normal  to  the  boundary  being  somewhat  larger 
than  the  range  of  the  elastic  field  of  the  boundary.  In  the 
unrelaxed  geometric  model  of  the  boundary  [Fig.  2(a)],  the 
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FIG.  I.  (a)  Construction  ofa  tilt-angle  grain  boundary  in  a  two  dimensional 
crystal.  Regions  A  and  6  show  the  upper  and  lower  halves  of  a  hexagonal 
bicrystal  Notice  the  arrays  in  one  region  are  rotated  by  an  angle  of  38.2 1*) 
with  respect  to  the  arrays  in  the  other  region.  Both  arrays  are  shown  in  C 
where  several  coincidence  sites  can  be  seen,  with  each  surrounded  by  a  ring 
of  atoms  from  either  A  or  B.  Throughout  this  paper,  the  kite-shaped  struc¬ 
ture  will  be  used  as  a  pattern  recognition  aid  outlining  the  arrrangement  of 
atoms  at  the  grain  boundary  The  boundary  itself  is  located  along  line  a-a 
Region  D  shows  a  portion  of  the  simulation  cell.  Ibl  Hexagonal  unit  cell  of 
DSC  lattice.  (c|  The  face-centered  rectangular  DSC  lattice  unit  cell  defining 
[100]  DSC  and  [010]  DSC 


density  of  the  boundary  region  is  the  same  as  that  in  a  perfect 
crystal.  There  is,  however,  a  pair  of  atoms  in  each  period, 
adjacent  to  the  kite  tail,  whose  spacing  is  appreciably  less 
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than  the  normal  nearest-neighbor  distance.  Such  a  “tight 
pair"  is  a  high  energy  configuration,  and  during  relaxation 
there  will  be  a  local  expansion  in  the  boundary  region  to 
relax  the  tight  pair.  This  relaxation  will  locally  decrease  the 
density  of  the  boundary  and  introduce  free  volume.  Unless 
the  cell  is  properly  adjusted  to  accommodate  the  free  vol¬ 
ume,  there  will  be  a  uniaxial  tensile  or  compressive  strain  in 
the  direction  normal  to  the  boundary,  superimposed  on  the 
usual  boundary  elastic  strain  field,  and  extending  through¬ 
out  the  cell. 

There  are  several  ways  of  adjusting  the  simulation  cell 
to  obtain  the  correct  free  volume.  We  have  tried  three  geo¬ 
metries,  all  of  which  yielded  essentially  the  same  results.  In 
all  cases,  we  monitored  the  potential  energy  and  spacing  of 
atoms  in  the  crystals  to  detect  the  presence  of  a  uniaxial 
strain.  In  two  calculations,  they  borders  were  periodic  and  in 
the  third  they  were  free  surfaces.  Although  we  finally  decid¬ 
ed  on  a  56-atom  cell  for  the  dynamic  simulation,  cells  with 
up  to  84  atoms  (YL  —  42  [010]  DSC)  were  investigated.  The 
most  definitive  results  were  obtained  on  the  84-atom  cell 
with  freey  borders,  but  the  results  of  the  other  calculations 
also  have  some  interesting  aspects. 

The  static  relaxations  were  carried  out  by  two  different 
simulation  techniques,  both  giving  essentially  the  same  re¬ 
sults.  In  the  molecular  dynamics  technique  the  atoms  moved 
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(a)  (b) 


FIG.  2.  Cell  of  56  atoms  used  in  the  present  computer  molecular  dynamics 
simulation  of  grain  boundary  migration,  (a)  Initial  unrelaxed  configuration 
showing  the  close  overlap  of  atoms  2  and  56,  and  29  and  30.  |b)  Configura¬ 
tion  prior  to  dynamic  simulation — fully  relaxed  including  relief  of  longitu¬ 
dinal  compression  by  slight  increase  of  cell  dimension  along  the  y  direction 
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under  the  influence  of  mutual  interactions,  but  as  the  system 
evolved  in  time,  a  small  damping  was  applied  so  that  the 
system  gradually  cooled  down  and  the  atoms  settled  into  a 
configuration  of  zero  net  force.13  Another  technique  used  to 
obtain  a  force-free  configuration  was  a  conjugate  gradient 
method.18 

In  the  calculations  with  periodic  y  borders  at  YL  =  28 
[010]  DSC,  the  56-atom  ceil  shown  in  Fig.  2(a)  was  relaxed 
without  any  changes  of  the  border  positions.  Upon  relaxa¬ 
tion  a  uniaxial  compressive  strain  developed.  In  order  to  re¬ 
move  this  strain,  we  increased  YL  in  small  steps  and  repeated 
the  relaxation  each  time.  This  procedure  led  to  a  minimum 
in  potential  energy  of  the  cell  at  YL  =  28.6085  [010]  DSC, 
and  an  average y  separation  of  0.5  [010]  DSC  for  the  ( 541 1 
planes  midway  between  the  boundary  (the  perfect  crystal 
value)  at  almost  this  same  value,  YL  =  28.6070  [010]  DSC. 
The  relaxed  structure  corresponding  to  the  latter  YL  value  is 
shown  in  Fig.  2(b).  This  structure  is  simply  the  initial  struc¬ 
ture  expanded  in  the  y  direction  by  0.304  [010]  DSC  per 
boundary,  i.e.,  giving  Vf  =  0.304  atomic  volume  per  period. 
Several  authors  have  stressed  the  importance  of  permitting 
global  rigid-body  translations  of  one  crystal  with  respect  to 
the  other  in  a  direction  parallel  to  the  boundary  (x  here)  as  a 
relaxation  adjustment.'9  The  present  relaxed  structure  for 
the  boundary  system  with  the  Lennard-Jones  potentials  did 
not  result  in  any  x  translation,  in  spite  of  the  fact  that  such 
translations  were  consistent  with  the  use  of  periodic  borders 
as  shown  by  the  calculations  to  be  described  next. 

In  a  second  calculation  using  periodic  y  borders,  we 
removed  one  of  the  atoms  in  the  tight  pair  in  each  boundary 
period.  This  created  a  54-atom  cell  with  YL  —  28  [010]  DSC, 
the  initial  length  of  the  56-atom  cell.  The  effect  of  a  series  of 
individual  atomic  relaxations  was  that  one  crystal  under¬ 
went  a  global  translation  in  the  x  direction,  with  respect  to 
the  other,  and  the  atoms  near  the  boundary  relaxed  in  to¬ 
ward  the  boundary  leaving  a  tensile  strain  in  the  center  of  the 
crystal.  Subsequent  adjustment  of  they  dimension  to  relieve 
this  tensile  strain  yielded  the  same  relaxed  configuration  as 
above  with  the  same  value  of  Vf  but  in  a  shorter  computa¬ 
tional  cell  corresponding  to  the  reduced  number  of  atoms. 
This  calculation  demonstrated  rather  dramatically  that  peri¬ 
odic  borders  do  not  inhibit  global  translations,  and  that  the 
structure  in  Fig.  2(a)  arrived  from  quite  different  paths  was 
very  likely  the  minimum  energy  structure.  It  also  shows  that 
the  structural  vacancy  arises  naturally  from  the  relaxation 
procedure  rather  than  from  an  arbitrary  initialization. 

In  both  the  54-  and  56-atom  cells,  adjustments  of  YL 
were  used  to  obtain  a  region  in  the  center  of  the  crystals 
which  on  the  average  had  no  uniaxial  strain.  However,  we 
could  not  eliminate  a  slight  strain  gradient  indicating  that 
the  elastic  fields  of  the  boundaries  were  overlapping  slightly. 
To  explore  this  we  increased  the  cell  size  to  84  atoms,  with 
YL  =  42  [010]  DSC,  and  used  free  surfaces  asy  borders.  This 
left  the  atoms  in  the  cell  quite  free  to  move  outward  so  that 
the  proper  free  volume  was  obtained  quite  naturally.  After 
the  relaxation  this  cell  had  a  region  7  [010]  DSC  in  length  in 
which  the  potential  energy  of  individual  atoms  was  within  1 
part  in  10*  of  the  values  obtained  in  the  other  calculations, 
and  the  interpianar  spacing  of  the  }  541 1  planes  was  within  1 
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FIG.  3.  (a)  The  initial  microstructure  produced  by  periodic  extension  in  the 
x  and  y  directions  of  the  computational  cell  ABCD.  (b)  The  microstructure 
after  one  boundary,  indicated  by  the  minus  sign,  has  migrated  downward. 


part  in  104  of  the  perfect  crystal  values.  To  the  same  degree  of 
resolution  the  strain  field  of  the  free  surface  extended  inward 
by  4.5  [010]  DSC,  and  that  of  the  boundary  outward  by  9.5 
[010]  DSC.  Vf  was  again  0.340  atomic  volume  period.  The 
extent  of  the  boundary  strain  field  in  units  of  the  boundary 
period  was  2.4.  This  is  comparable  to  what  we  found  in  relax- 
ingaX  =  5  [310]  twin  boundary  in  three-dimensional  simu¬ 
lations  of  a  boundary  in  a  fee  crystal  with  Lennard- Jones 
potential.20 

It  might  naively  have  been  expected  that  a  grain  bound¬ 
ary  would  perturb  the  crystal  away  from  it  over  a  smaller 
distance  than  would  a  free  surface.  The  present  results  to  the 
contrary  suggest  a  perhaps  equally  naive  explanation  that 
because  the  free  surface  is  unconstrained,  local  relaxation 
can  take  place  to  a  greater  degree  than  in  the  grain  boundary 
without  requiring  as  strong  a  perturbation  below  the  surface. 

The  above  results  show  that  in  the  56-atom  cell  used  in 
the  molecular  dynamics  simulations  the  boundaries  perturb 
each  other  through  their  elastic  fields.  Initially  the  micro- 
structure  generated  by  the  periodic  extension  of  the  cell  in 
the  x  and  y  directions  consists  of  an  array  of  equally  spaced 
boundaries  [Fig.  3(a)].  When  the  symmetry  is  broken  by 
thermal  fluctuations,  a  pair  of  boundaries  may  migrate  to¬ 
gether  with  continuous  attraction,  continually  “consuming” 
the  crystalline  region  in  between,  until  they  finally  annihilate 
each  other  leaving  in  the  system  a  single  crystal.  In  the  mo¬ 
lecular  dynamics  simulation  on  the  56-atom  cell,  this  pro¬ 
cess  was  actually  observed  at  sufficiently  high  temperatures. 

It  should  be  noted  that  as  a  final  check  on  our  static 
relaxation,  we  ran  the  molecular  dynamics  program  without 
any  damping  force  and  assigned  all  the  particles  zero  initial 
velocity.  If  the  system  were  not  in  a  relaxed  configuration  the 


particles  would  experience  a  net  force  and  therefore  would 
move.  With  the  configuration  shown  in  Fig.  2(b)  the  particles 
remained  essentially  motionless  for  several  hundred  time 
steps  of  simulation. 

III.  THRESHOLD  FOR  GRAIN-BOUNDARY  SLIDING 

In  any  consideration  of  the  dynamics  of  grain  boundar¬ 
ies,  sliding  is  one  of  the  basic  modes  of  motion.21-23  When 
viewed  in  terms  of  atomic  displacements  this  is  a  highly  co¬ 
operative  process  since  it  involves  the  translation  of  one 
group  of  particles  relative  to  another  group.  To  identify  such 
processes  in  our  simulation  we  have  adopted  the  following 
convention  which  seems  reasonable  to  us  but  is  admittedly 
arbitrary.  Since  instantaneous  atomic  displacements  have  a 
component  associated  with  thermal  vibrations  present  in 
any  finite-temperature  simulation,  we  consider  only  dis¬ 
placements  that  have  been  averaged  over  100  time  steps,  an 
interval  that  is  several  times  the  typical  vibrational  periods. 
With  this  averaging  we  expect  that  much  of  the  short-time 
vibrational  effects  are  eliminated.  Futhermore,  we  will  re¬ 
gard  sliding  to  have  occurred  when  the  relative  displacement 
of  one  group  of  atoms  relative  to  another  is  equal  to  or 
greater  than  one  DSC  unit  along  the  x  direction. 

Figure  4  shows  a  typical  pattern  of  average  particle  dis¬ 
placements  in  the  x  direction  during  a  relatively  early  stage 
of  a  simulation.  One  can  see  that  particles  numbered  1 
through  25  all  have  positive  displacements.  This  configura¬ 
tion,  which  indicates  that  boundary  sliding  has  occurred, 
developed  after  the  system  was  allowed  to  come  to  equilibri¬ 
um  at  a  particular  temperature.  Boundary  sliding,  however, 
did  not  always  occur  whatever  the  temperature;  in  fact,  for  a 
given  length  of  simulation  there  was  a  temperature  below 
which  sliding  was  not  observed. 

Since  all  simulations  began  with  a  relaxed  configura¬ 
tion,  the  energy  necessary  to  activate  the  sliding  process 
must  come  from  the  thermal  energy  of  the  particles.  While  it 
is  difficult  to  define  and  determine  precisely  the  potential 
energy  barrier  to  sliding,  we  have  established  a  correlation 
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FIG.  4.  Atomic  displacements  along  the  x  direction  |in  units  of  XDSC).  The 
displacements  are  averaged  values  over  an  interval  of  100  time  steps  and  are 
taken  during  a  relatively  early  stage  of  the  simulation. 
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FIG.  5.  Variation  of  “waiting  time”  (in  units  of  simulation  time  steps}  prior 
to  the  onset  of  boundary  sliding  as  a  function  of  the  average  kinetic  energy 
per  particle  (in  unit  of  e,  the  well  depth  of  the  Lennard-Jones  potential 
function). 


between  the  average  kinetic  energy  per  particle  and  the  time 
at  which  sliding  commences.  Figure  5  shows  that  at  lower 
system  temperatures  the  waiting  time  for  sliding  increases. 
This  behavior  is  quite  reasonable  since  the  activation  energy 
is  provided  by  thermal  fluctuations  and  the  necessary  fluctu¬ 
ations  should  occur  less  frequently  at  lower  temperatures. 
As  also  shown  in  Fig.  5,  at  sufficiently  high  temperatures  the 
time  at  which  sliding  occurs  appears  to  be  insensitive  to  tem¬ 
perature,  thus  indicating  that  even  when  there  is  sufficient 
thermal  energy  to  overcome  the  barrier  the  system  requires  a 
certain  period  of  time  for  equilibration,  perhaps  dependent 
on  dynamical  processes.  At  an  average  kinetic  energy  per 
particle  of  0.02 5*  sliding  was  not  observed  in  a  simulation 
run  extending  out  to  IS  000  time  steps.  We  regard  this  as 
evidence  that  the  sliding  processes  observed  here  are  ther¬ 
mally  activated.  In  such  thermally  activated  processes  a 
threshold  temperature  exists  below  which  the  probability  of 
sliding  within  a  given  time  interval  decreases  exponentially. 

IV.  COUPLED  BOUNDARY  SLIDING  AND  MIGRATION 

A  competing  process  to  boundary  sliding  is  boundary 
migration. 2,  22  This  process  involves  the  motion  of  the  grain 
boundary  in  they  direction  in  our  simulation  cell.  Our  simu¬ 
lation  results  showed  that  both  migration  and  sliding  pro¬ 
cesses  could  take  place,  and  the  general  movements  of  a 
grain  boundary  can  be  described  as  sequences  of  coupled 
sliding  and  migration.  It  is  possible  to  follow  the  detailed 
•  sequence  of  atomic  motions  which  occur  during  migration. 
These  motions  involve  the  correlated  sequential  motion  of  a 
number  of  atoms.  Two  different  sequences  have  been  ob- 
’  served  in  which  the  boundary  moves  either  three  or  four 
{ 341 )  planes  during  the  migration  sequence.  As  will  be  dis¬ 
cussed  in  more  detail  later,  the  migration  sequences  also  in¬ 
volve  a  localized  shearing  or  sliding  of  one  crystal  with  re¬ 
spect  to  the  other.  The  sequences  are  in  fact  those  of 
migration-sliding. 

By  following  a  group  of  atoms  (two  groups  are  indicated 
in  heavy  circles),  one  can  readily  discern  how  the  structure 
changes  in  the  process.  In  the  initial  configuration,  Fig.  6(a), 
a  kite-shaped  group  shows  the  boundary  plane  is  at  atom  1 


and  a  ring  of  atoms  surrounds  atom  3  which  is  part  of  the 
upper  bicrystal.  At  the  end  of  the  sequence,  Fig.  6(c),  the 
boundary  has  moved  to  atom  4,  a  migration  of  three  { 341 ) 
planes,  while  atom  3  is  now  clearly  part  of  lower  bicrystal. 
Figure  6(b)  shows  a  saddlepoint  or  activated  state  configura¬ 
tion  with  atom  3  in  the  middle  of  the  ring. 
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FIG.  6.  Sequence  of  upper  (rain  boundary  configuration  as  in  Fig.  2(b) 
showing  coupled  sliding  and  migration.  Each  configuration,  except  that  of 
r  -  0,  is  an  average  over  an  interval  of  100  time  steps  starting  at  time  step.  >. 
(a)  t  -  0,  (b)  f  =  1401,  (c|r=  1501,(d|f  «  1701. 
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Comparing  the  vertical  alignment  of  atom  pain  (2,5) 
and  (4,1)  at  the  beginning  and  the  end  of  the  sequence,  one 
can  detect  a  sliding,  or  shearing,  of  the  upper  bicrystal  to  the 
right.  Thus  in  the  same  sequence  sliding  and  migration  have 
both  occurred.  Notice  that  atom  6  was  lined  up  initially  with 
atom  52,  and  at  the  end  of  the  sequence  it  was  lined  up  with 
atom  55.  The  amount  of  sliding  is  [100]  DSC  by  referring  to 
the  positions  of  atom  52  and  55  in  Fig.  1. 

We  have  also  observed  and  carried  out  a  similar  analy¬ 
sis  of  migration-sliding  sequence  involving  the  boundary  tra¬ 
versing  four  j  541 1  layers  (2  [010]  DSC).  The  activated  state 
is  the  same,  but  the  relative  shear  is  different. 

V.  GRAIN-BOUNDARY  ANNIHILATION 

The  two  grain  boundaries  in  the  simulation  cell  are  cou¬ 
pled  through  the  strain  field  and  through  the  thermal  vibra¬ 
tions.  If  the  boundaries  move  toward  each  other  sufficiently, 
they  can  combine  to  form  a  perfect  crystal.  We  will  call  this 
process  “boundary  annihilation.”  One  can  calculate  the  po¬ 
tential  energy  of  the  system  with  two  boundaries  in  its  re¬ 
laxed  configuration,  as  shown  in  Fig.  2(b),  and  the  potential 
energy  of  a  perfect  crystal  specimen  with  the  same  number  of 
particles  and  the  same  area.  This  difference  in  potential  ener¬ 
gy  contributes  to  the  activation  energy  for  boundary  annihi¬ 
lation.  However,  boundary  annihilation  will  not  occur  spon¬ 
taneously  at  zero  temperature  since  sliding  and  migration 
processes  require  thermal  activations. 

We  have  observed  boundary  annihilation  in  our  simula¬ 
tions.  The  annihilation  process  is  apparently  quite  compli¬ 
cated,  and  it  did  not  occur  immediately  after  the  commence¬ 
ment  of  sliding  and  migration.  The  picture  which  emerged 
from  our  results  is  that  at  a  given  temperature  above  the 
threshold  the  system  requires  a  certain  average  time  interval, 
which  we  may  call  “incubation  period,”  to  initiate  the  slid¬ 
ing  and  migration  processes.  Then  there  follows  a  time  inter¬ 
val  during  which  the  boundaries  move  in  a  stochastic  man¬ 
ner.  Finally  the  two  grain  boundaries  steadily  move  toward 
each  other  and  annihilate. 

The  post  annihilation,  single  crystal  structure  is  shown 
in  Fig.  7(c).  From  this  result  and  the  initial  configuration, 
Fig.  7(a),  with  two  grain  boundaries  one  obtains  a  typical 
displacement  field  giving  the  net  atomic  displacements.  Sub¬ 
tracting  the  atomic  positions  in  the  initial  configuration  one 
obtains  the  indicated  net  displacement  field.  From  this,  one 
can  infer  that  the  boundaries  did  not  migrate  through  the 
central  region  of  the  simulation  cell  on  their  way  to  annihila¬ 
tion.  Details  of  the  displacement  field  clearly  depend  on  the 
particular  sequence  of  sliding  and  migration  processes  in¬ 
volved,  and  different  field  patterns  have  been  obtained  in 
simulation  runs  using  different  sets  of  random  numbers  for 
initialization. 

It  is  instructive  to  examine  the  annihilation  process 
from  the  standpoint  of  the  excess  thermodynamic  functions 
per  unit  area  of  the  grain  boundary.  The  excess  energy  arises 
from  the  different  average  configuration  of  the  boundary 
atoms  with  respect  to  those  in  the  crystal,  as  well  as  the 
differences  between  the  contribution  to  potential  energy 
arising  from  the  temperature-dependent  position  fluctu¬ 
ations.  This  latter  contribution  results  from  the  different  an- 
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FIG.  7.  Simulation  cell  showing  the  actual  net  motions  of  the  atoms  asso¬ 
ciated  with  the  migration  of  the  two  grain  boundaries  from  initial  configura¬ 
tion  {a)  to  annihilation.  The  total  net  displacements  (magnified  by  a  factor  of 
2)  of  the  atoms  are  given  in  (b).  while  (c)  represents  the  atomic  configuration 
averaged  over  100  time  steps  after  annihilation. 


harmonic  vibrational  energies  of  the  boundary  atoms  from 
those  of  atoms  in  the  perfect  lattice.  There  is  also  an  excess 
entropy  arising  from  corresponding  differences  in  harmonic 
and  anharmonic  vibrational  entropies  between  the  grain¬ 
boundary  and  the  perfect-lattice  systems.  Finally,  there  may 
be  a  PA  V  contribution  to  the  excess  enthalpy  or  Gibbs  free 
energy  arising  from  the  excess  boundary  volume.  The  anni¬ 
hilation  process  was  observed  to  produce  an  increase  of  aver¬ 
age  temperature,  or  equivalently  the  average  kinetic  energy 
per  particle,  as  shown  in  Fig.  8.  Annihilation  occurred  in  this 
case  at  about  time  step  3600. 

VI.  EXCESS  FREE  AREA  AND  STRAIN  DISTRIBUTION 

We  have  thus  far  considered  the  dynamics  of  grain 
boundaries  in  terms  of  individual  particle  motions.  Another 
way  of  discussing  the  simulation  results  is  to  analyze  the 
distribution  of  excess free  area  in  the  system  at  various  stages 
of  dynamical  evolution.  Our  simulation  cell  containsa  total 
excess  free  area  equal  to  1.22  \  V/  =  0.304,  two  boundaries 
and  two  periods)  times  the  area  occupied  by  an  atom  in  the 
perfect  crystal  when  we  compare  the  cell  area  of  the  bound¬ 
ary  system  to  that  of  the  corresponding  perfect  crystal.  Most 
of  this  excess  is  initially  localized  in  the  immediate  region  of 
the  two  grain  boundaries,  but  as  the  boundaries  begin  to  slide 
and  migrate,  the  free-area  distribution  can  also  change. 

To  determine  the  distribution  of  free  area  under  dy¬ 
namical  conditions,  we  have  adopted  the  procedure  of  Ar¬ 
gon  and  Kuo24  for  estimating  the  excess  free  area  in  a  local 
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FIG.  8.  Instantaneous  and  time-averaged 
kinetic  energy  per  particle  (in  unit  of  poten¬ 
tial  well  depth  f).  One  can  observe  an  in¬ 
crease  in  kinetic  energy  after  time  step  3600 
as  a  result  of  boundary  annihilation.  Also 
shown  is  the  root  mean-square  deviation. 


region  (Fig.  9).  While  the  excess  free  area  thus  defined  is  not 
always  entirely  concentrated  in  the  grain-boundary  regions, 
movements  of  the  excess  free  area  indeed  follow  the  motions 
of  the  boundaries.23  When  the  boundaries  annihilate  at 
approximately  time  step  3600,  the  excess  free  area  cannot 
vanish  because  the  simulation  is  carried  out  at  constant  area. 
Moreover,  it  appears  that  the  excess  free  area  is  not  redistri¬ 
buted  uniformly  throughout  the  cell;  rather,  it  remains  local¬ 
ized  though  not  static.  It  is  possible  that  the  free  areas  can 
coalesce  and  form  a  vacancy,  although  we  have  not  observed 
this  process  in  any  of  our  simulations. 

In  an  attempt  to  display  the  deformations  that  occur 
during  sliding  and  migration,  we  have  carried  out  an  analysis 
of  the  local  shear  and  dilation  strains  using  a  procedure24 
that  makes  use  of  particle  positions  at  two  successive  time 
steps.  The  results  are  shown  in  Fig.  10,  where  the  strains 
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FIG  9.  DiM  rt  but  ion  of  excess  free  area  in  ihc  simulation  cell  under  dynamic 
conditions.  The  cell  is  divided  into  540  boxes  with  O0  the  y  direction 
At  a  fixed  y  coordinate,  only  the  box  with  the  largest  free  area  is  considered, 
and  a  symbol  is  shown  for  that  v  coordinate  if  this  free  area  exceeds  the 
average  free  area  value  in  the  corresponding  perfect  crystal  cell  Three  sym¬ 
bols  are  used  to  indicate  relative  magnitudes  *,  and  l  in  increasing 
order  of  excess  free  area. 
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FIG.  10.  Strain  distribution  in  the  simulation  cell  under  dynamic  condi¬ 
tions  The  cell  is  divided  into  540  boxes  with  60  along  the  y  direction  Ai  a 
fixed  y  coordinate,  only  the  box  with  the  largest  strain  is  considered,  and  a 
symbol  is  shown  for  the  >■  coordinate  if  the  strain  exceeds  the  average  strain 
in  the  contspnnding  perfect  crystal  cell.  Three  symbols  are  used  to  indicate 
relative  magnitudes  « and  1.  in  decreasing  values  of  strain  la|  Dilation 
strain,  (b|  shear  strain. 
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have  been  averaged  over  100  time  steps  and  different  sym¬ 
bols  are  used  to  indicate  relative  magnitudes  as  was  done  in 
Fig.  9.  Both  shear  and  dilation  strains  are  seen  to  be  higher  in 
the  grain-boundary  regions,  as  we  would  expect.  When  there 
was  no  boundary  migration,  their  magnitudes  were  relative¬ 
ly  small.  During  migration  high  strains  were  developed  in 
the  boundary  region,  indicating  that  when  atoms  have  to 
slide  relative  to  each  other  they  also  induce  considerable  de¬ 
formations  in  the  perpendicular  direction. 

In  Fig.  10  one  can  see  that  although  the  system  was 
initially  strain  free,  both  shear  and  dilation  strains  quickly 
built  up  after  the  start  of  the  simulation,  and  they  moved 
through  the  cell  following  the  grain  boundaries.  Also,  the 
strains  were  not  concentrated  only  in  the  boundary  regions. 
For  example,  around  time  steps  2500-3000,  there  appeared  a 
region  of  shear  and  dilation  strains  in  the  center  of  the  cell 
while  the  two  boundaries  were  actually  located  at  the  top 
and  bottom  edges  of  the  cell.  After  annihilation  there  was 
considerable  residual  strain  in  both  shear  and  dilation. 
These,  however,  were  not  distributed  in  any  localized  man¬ 
ner. 
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-'Our  simulation  cell  is  such  that  if  we  divide  the  cell  area  into  60  squares.  3 
along  the  x  direction  and  20  along  they  direction,  there  will  be  roughly 
one  particle  per  square.  To  obtain  some  spatial  resolution  we  then  divide 
the  cell  into  540  boxes,  9  along  thex  direction  and  60along  they  direction 
In  addition,  each  box  is  further  divided  into  25  small  boxes.  5  along  each 
direction.  At  any  given  time  the  small  boxes  which  contain  the  atom 
centers  are  found,  these  boxes  and  the  224  nearest  boxes  surrounding  each 
one  are  considered  to  be  occupied  (overlap  is  allowed).  By  counting  the 
number  of  unoccupied  small  boxes  in  each  large  box  we  determine  the  free 
area  in  that  large  box.  We  then  repeat  the  same  procedure  for  the  perfect 
crystal.  By  excess  free  are  we  mean  those  cases  where  the  number  of  unoc¬ 
cupied  small  boxes  exceeds  the  maximum  number  of  unoccupied  small 
boxes  found  in  the  perfect  crystal.  These  cases  are  then  divided  into  three 
groups  with  each  indicated  by  a  different  symbol  in  Fig.  12 
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